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Abstract 

The  plotting  of  streamlines  is  an  effective  way  of  visualizing  fluid  motion  in  steady  flows.  Ad¬ 
ditional  information  about  the  flowfield,  such  as  local  rotation  and  expansion,  can  be  shown 
by  drawing  in  the  form  of  a  ribbon  or  tube.  In  this  paper,  we  present  efficient  algorithms  for 
the  construction  of  streamlines,  streamribbons  and  streamtubes  on  unstructured  grids.  A  spe¬ 
cialized  version  of  the  Runge-Kutta  method  has  been  developed  to  speed  up  the  integration  of 
particle  paths.  We  have  also  derived  closed-form  solutions  for  calculating  angular  rotation  rate 
and  radius  to  construct  streamribbons  and  streamtubes,  respectively.  According  to  our  analysis 
and  test  results,  these  formulations  are  two  to  four  times  better  in  performance  than  previous 

numerical  methods.  As  a  large  number  of  traces  are  calculated,  the  improved  performance  could  _ 
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1  Introduction 


Streamlines,  streamribbons  and  streamtubes  are  very  powerful  techniques  for  visualizing  steady 
vector  fields.  A  streamline  is  the  path  of  a  massless  particle  which  is  released  in  a  steady  flow. 
The  plotting  of  the  particle  paths  produces  a  streamline  picture,  which  is  of  both  qualitative 
and  quantitative  value  to  the  engineer.  Streamline  pictures  allow  the  engineer  to  visualize  fluid 
motion  and  to  locate  regions  of  high  and  low  velocity  and,  from  these,  zones  of  high  and  low 
pressure. 

Given  a  fluid  flow  with  velocity  field  u(x(t )),  a  streamline  is  an  integral  curve  of  u.  That  is,  a 
streamline  can  be  calculated  by  solving  the  following  equation: 

AS  =  5(S(«))  (1) 

where  t  is  a  parameter  along  the  streamline  and  is  not  to  be  confused  with  time  [11]. 

A  streamribbon  can  show  the  translation,  angular  rotation,  and  rates  of  shear  deformation  of 
the  flow.  Ideally,  it  is  constructed  by  tracing  a  set  of  streamlines  originated  from  multiple  seed 
locations  on  a  straight  line  segment.  That  is,  the  path  swept  by  the  deformable  line  segment 
becomes  a  streamribbon.  Volpe  [12]  constructs  a  streamribbon  in  this  fashion  by  tracing  a  large 
number  of  adjacent  streamlines.  However,  the  number  of  streamlines  needed  to  form  smooth 
ribbon  surfaces  could  be  tremendous  and  the  corresponding  computational  cost  would  be  high. 
In  practice,  the  construction  of  streamribbons  is  simplified,  though  some  information  such  as 
shear  deformation  would  be  lost.  In  [4],  a  streamribbon  is  generated  by  computing  only  a 
few  streamlines  and  creating  polygons  between  adjacent  streamlines  to  form  the  surface  of  the 
streamribbon.  This  method  still  requires  complicated  algorithms  to  deal  with  the  convergence, 
the  divergence  and  the  splitting  of  streamribbons.  Darmofal  and  Haimes  [2],  Ma  and  Smith 
[7],  and  Pargendarm  [9]  use  one  streamline  and  vectors  normal  to  the  local  velocity  to  form  a 
streamribbon.  In  this  way,  the  resulting  ribbons  only  show  the  translation  and  angular  rotation 
of  the  flow.  We  adopt  Darmofal  and  Haimes’  algorithm  by  using  two  parallel  edges  to  form  a 
streamribbon.  First,  a  streamline  is  generated  to  serve  as  the  first  edge  of  the  streamribbon.  A 
normal  vector  is  calculated  at  each  point  of  the  streamline  by  rotating  a  constant  length  vector 
about  the  streamline.  Then  the  second  edge  of  the  streamribbon  is  formed  by  connecting  the 
end  points  of  the  normal  vectors. 

Formally,  a  streamtube  is  defined  as  the  surface  formed  by  all  streamlines  passing  through 
a  given  closed  curve  in  the  flow  [11].  Streamtubes  are  used  to  visualize  expansion,  contraction 
and  deformation  of  the  flow.  In  [2],  a  streamtube  is  created  by  connecting  the  circular  crossflow 
sections  along  a  streamline.  The  radius  of  a  cross  flow  section  is  determined  by  the  local  cross 
flow  expansion  rate.  A  streamtube  constructed  in  this  manner  does  not  reveal  the  deformation 
of  the  flow.  Again,  this  is  a  technique  more  computationally  feasible  and  we  adopt  it  in  this 
work.  In  [7],  to  visualize  both  flow  convection  and  diffusion,  statistical  dispersion  of  the  fluid 
elements  about  a  streamline  is  computed  by  using  added  scalar  information  about  the  root  mean 


square  value  for  the  vector  field  and  its  Larangian  time  scale.  The  result  defines  the  radius  of  the 
cross  flow  section  and  also  forms  a  tube-like  surface.  Schroeder  et.  al.  [10]  introduce  a  technique 
called  Stream  Polygon  for  visualizing  local  deformation  of  the  flow. 

In  this  paper,  we  present  efficient  algorithms  to  compute  streamlines,  streamribbons  and 
streamtubes  on  unstructured  grids.  Our  algorithms  are  mainly  based  on  those  developed  in 
[2].  Several  new  computational  techniques  are  derived  and  used  to  improve  performance.  These 
new  computational  techniques  include  a  specialized  version  of  Runge-Kutta  method,  a  simpler 
procedure  to  compute  the  angular  rotation  rate  of  the  flow  and  an  explicit  solution  for  calcu¬ 
lating  the  radius  of  streamtube.  The  overview  of  our  algorithms  is  described  in  Section  2.  The 
new  computational  techniques  are  derived  in  Section  3.  Data  structures  used  and  the  memory 
requirements  for  using  the  algorithms  are  described  in  Section  4.  Finally,  we  present  some  ex¬ 
perimental  results  using  three  different  data  sets  to  demonstrate  the  time  efficiency  of  the  new 
particle  tracing  algorithm. 


In  this  paper,  we  assume  that  all  cells  are  linear  tetrahedra.  Other  types  of  cells  have  to  be 
decomposed  into  tetrahedra  in  preprocessing  stages.  In  a  tetrahedral  cell,  the  three  components 
of  the  vector  field  are  linear  functions  of  the  physical  coordinates.  Their  interpolation  functions 
can  be  formulated  as: 

ui(x,y,z)  =  (tiX  +  biy  +  cx  z  +  d1, 

u2(x,y,z)  =  a2x +  b2y  +  c2z  +  d2,  (2) 

u3(x,  y ,  z)  —  a3x  +  b3y  +  c3z  +  d3. 

where  u;,  z  =  1,2,3,  are  the  three  components  of  vector  field;  at,  b *,  ct,  di,i  —  1, 2, 3,  are  the  coef¬ 
ficients  of  the  interpolation  functions;  x,y,z  are  the  physical  coordinates.  The  above  equations 
can  be  re-written  in  a  concise  form: 
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When  calculating  a  streamline,  it  is  necessary  to  find  the  cell  in  which  this  streamline  enters 
at  each  time  step.  A  method  is  given  in  [6]  to  solve  this  problem.  In  this  method,  the  physical 
coordinates  of  the  point  calculated  at  each  time  step  are  transformed  into  the  canonical  coor¬ 
dinates  as  shown  in  Figure  1.  Then  the  canonical  coordinates  are  used  to  determine  the  cell 
which  the  streamline  enters.  Since  the  cells  are  linear  tetrahedra,  the  transformation  functions 
of  the  canonical  coordinates  are  linear  functions  of  the  physical  coordinates.  That  is,  we  use  the 
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Transformation 


Physical  Coordinates  Canonical  Coordinates 

Figure  1:  Coordinate  System  Transformation 

following  formulation  to  convert  the  physical  coordinates  into  the  canonical  coordinates: 

£  =  Rx  +  k 


where  x  is  a  physical  coordinate  vector,  f  is  the  canonical  coordinate  vector  of  x .  and 
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To  obtain  r/j  and  ki,  we  can  solve  the  following  linear  system: 
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where  aq,  yi  and  are  the  physical  coordinates  of  each  node  of  the  the  cell. 

2.1  Streamline  Construction 

Given  an  initial  point  in  a  physical  domain,  a  streamline  can  be  calculated  by  solving  Equa¬ 
tion  1.  The  4th  order  Runge-Kutta  method  is  applied  to  integrate  the  equation  stepwise.  After 
calculating  a  point  of  the  streamline,  Equation  7  is  used  to  transform  the  physical  coordinates 
of  the  point  into  the  canonical  coordinates.  If  all  the  three  components  of  the  canonical  coordi¬ 
nates  are  between  0.0  and  1.0,  this  point  is  still  inside  the  current  cell  where  the  computation  of 
the  point  takes  place.  The  coefficients  of  the  interpolation  functions  of  the  current  cell  are  still 
valid  for  next  step  integration.  Otherwise,  a  searching  for  a  new  cell  which  contains  the  point 
is  started  according  to  the  canonical  coordinates.  After  finding  the  new  cell,  the  computation 
of  next  position  can  be  performed.  This  pattern  of  calculation  is  repeated  until  the  streamline 
reaches  a  physical  boundary  or  the  number  of  time  steps  exceeds  a  pre-defined  limit. 
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Normal  Vector 


Figure  2:  Example  of  Streamribbon  Construction 


A  streamribbon  has  two  edges  as  we  have  described.  The  first  edge  of  a  streamribbon  is  con¬ 
structed  by  calculating  a  streamline,  and  the  second  edge  is  generated  by  connecting  the  end 
points  of  the  normal  vectors  of  the  streamline.  The  normal  vectors  are  calculated  by  rotating 
a  constant  length  vector  about  the  streamline  at  each  point  of  the  streamline.  The  constant 
length  vector  can  be  any  vector  which  is  orthogonal  to  the  streamline  at  the  initial  point.  The 
surface  of  the  streamribbon  is  then  formed  by  connecting  the  end  points  of  the  normal  vectors 
and  their  corresponding  points  on  the  streamline.  An  example  is  depicted  in  Figure  2.  The  angle 
of  rotating  the  constant  length  vector  is  governed  by: 


d0_ 

dt 

w 


curliu) 


(11) 

(12) 


where  9  is  the  rotation  angle.  Equations  11  and  1  are  solved  stepwise  when  constructing  a 


streamribbon. 


A  streamtube  is  created  by  generating  a  streamline  and  by  connecting  the  circular  crossflow 
sections  along  the  streamline.  The  radius  of  a  streamtube  is  governed  by  the  following  ordinary 
differential  equation: 


1  dr 
r  dt 

Vx  •  u 


(14) 

(15) 


where  r  is  the  streamtube  radius,  Vx  •  u  is  the  local  cross  flow  divergence,  and  ^7  represents  the 
change  of  velocity  magnitude  along  the  streamline.  Equations  1,11  and  14  are  solved  stepwise 
when  constructing  a  streamtube.  Equation  1  is  used  to  calculate  the  center  of  the  streamtube, 
while  Equations  11  and  14  are  used  to  calculate  the  angle  of  rotation  and  the  radius  of  the 
streamtube.  Figure  3  contains  an  example  of  constructing  a  streamtube. 
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Figure  3:  Example  of  Streamtube  Construction 


3  New  Computational  Methods 

In  order  to  construct  streamlines,  streamribbons  and  streamtubes,  we  need  to  solve  the  ODE’s 
mentioned  in  the  previous  sections.  Based  on  the  interpolation  functions  of  linear  tetrahedral 
cells,  we  have  developed  specialized  ODE  solvers  to  speed  up  our  algorithms. 

3.1  A  Specialized  Version  of  the  Runge-Kutta  Method  for  Stream¬ 
line  Construction 

By  combining  Equations  1  and  3  the  governing  equation  of  a  streamline  can  be  formulated  as: 


=  f(x,  t)  =  Bx  +  d  (16) 

The  4th  order  Runge-Kutta  method  is  applied  to  solve  this  ODE: 

x(t  +  h)  =  x{t)  +  ^(F1  +  2F2  +  2F3  +  F4)  (17) 

F4  =  hf(x,t )  (18) 

F2  =  hf(x  +  F1/2,t  +  h/2)  (19) 

F3  =  hf(x  +  F2/2,t  +  h/2)  (20) 

F4  =  hf(x  +  F3,t  +  h)  (21) 


where  h  is  the  time  step  size.  By  substituting  Equations  16  and  3  into  the  right  hand  sides, 
Equations  18  -  21  can  be  expanded  as: 

Fi  =  hf(x ,  t) 

—  h(Bx  +  d ) 

F2  =  hf{x  +  Ei/2,  t  +  h/2) 

=  h(B(x  +  F1/2)  +  d) 

=  (h2B/2  +  h){Bx  +  d) 

F3  —  hf(x-)-F2/2,t-\-h/2) 

=  h(B(x  +  F2/2')  +  d) 


5 


=  (h3B2/4:  +  h2B/2  +  h){Bx  +  d) 

F4  =  hf(x  +  F3,t  +  h) 

=  h(B(x  +  E3)  +  d) 

=  (h4B3/4:  +  hzB2/ 2  +  h2B  +  h)(Bx  +  d) 


By  using  these  equations,  the  Runge-Kutta  method  shown  in  Equation  17  can  be  expressed  as: 


x(t  +  h)  =  x(t)  +  ~(jFi  +  2E2  +  2E3  +  F4 ) 


/r  (AB)2  (d£)3  (hB)4.  t>. 

(I  +  hB  + 

lh(r  .  hB  \hBf  \hBf  - 

+i(/+ir  +  _i_  +  -4rM 


=  +  tf2d  (22) 

=  #!*(<)  +  £  (23) 

Since  B  and  d  are  constants,  H\  and  d\  can  be  calculated  by  using  Horner’s  algorithm  [3]. 
Hereafter  the  computations  of  the  4th  order  Runge-Kutta  method  require  only  a  matrix-vector 
multiplication  and  a  vector-vector  addition. 


3.2  Explicit  Solution  for  the  Angular  Rotation  Rate 

The  angular  rotation  rate  is  governed  by  the  ODE  formulated  in  Equation  11.  Since  the  velocity 
u  is  linear  within  a  cell,  the  curl  of  u  is  a  constant  vector.  According  to  Equation  3,  we  have: 


w  =  curl(u) 

-  V  x  u 

=  V  x  ( Bx  +  d) 

=  [  b3  -  c2  ci  — az  a2 
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where  0(0)  is  the  rotation  angle  at  the  previous  time  step,  6(h)  is  the  rotation  angle  at  the  current 
time  step,  u(h)  is  the  velocity  at  the  current  time  step,  and  u(0)  is  the  velocity  at  the  previous 
time  step.  This  closed  form  solution  is  used  to  compute  the  rotation  angle  of  the  normal  vector 
about  the  streamline.  The  only  unknown  values  involved  in  this  solution  are  u(h )  and  its  velocity 
magnitude.  Since  u(h )  can  be  calculated  by  using  Equation  3,  the  major  cost  of  this  solution  is 
reduced  to  a  matrix- vector  multiplication. 


3.3  Explicit  Solution  for  the  Radius  of  Streamtube 

The  governing  equation  of  streamtube  radius  is  shown  in  Equation  14.  This  ODE  can  be  solved 
analytically: 

rh  1  .  1  fh  „ 

/  -dr  =  -  Vx  •  udt 

Jo  r  2  Jo 

In  (rh.)  -  ln(r0)  =  \[  Vx  ■  udt 

i  Jo 

Mr,)  =  Mr„)  +  i (J‘  v-aat-J"  ^ it) 

From  Equations  3  and  1,  the  divergence  of  u  is: 

V  ■  u  =  di  +  62  T  C3 

and 


dx  =  u  dt 


Therefore, 


In(r^)  = 
rh  = 


1  (I'll 

ln(r0)  +  -((a!  +  b2  +  c3)h  -  /  — r) 
z  Jo  u 

1 

r0exp(-(ai  +  b2  +  c3)h  -  \n(uh)  +  ln(tt0)) 


(24) 

(25) 


i  'LL 

rh  =  roexp(-(aa  +  b2  +  c3)h)-y- 


Ul 


(26) 


Equation  26  is  used  to  compute  the  radius  of  streamtube,  where  is  the  streamtube  radius  at 
the  current  time  step,  ro  is  the  radius  at  the  previous  time  step,  u0  and  uh  are  the  magnitudes 
of  velocity  at  the  previous  step  and  the  current  step.  Since  the  magnitude  of  velocity  at  current 
step  has  been  calculated  when  computing  the  angle  of  rotation,  there  is  no  unknown  value  in  the 
right  hand  side  of  this  equation.  The  cost  of  calculating  rh  composes  only  a  few  multiplications. 


3.4  Integration  Step  Size 

The  value  of  h  is  crucial  for  integration  of  particle  paths.  In  [1],  Buning  suggested  to  choose 
this  time  step  size  based  on  the  cell  size  and  the  inverse  of  velocity  magnitude.  Darmofal  [2] 
used  a  similar  method  to  determine  the  value  of  h  for  tracing  particle  paths,  but  for  constructing 
streamribbons,  h  is  further  restricted  by  the  angle  of  rotation  to  produce  a  smoother  ribbon 
surface.  In  our  current  implementation,  h  is  fixed  for  the  entire  streamline.  A  default  step  size  is 
determined  for  the  overall  domain  by  using  Buning’s  method  at  the  preprocessing  stage,  though 
h  can  be  interactively  modified. 

4  Data  Structures 

To  implement  the  above  methods,  the  major  data  structures  are  composed  of  a  list  of  cell  records 
and  a  list  of  node  records.  To  further  speed  up  the  construction  of  streamlines,  streamribbons 
and  streamtubes,  at  the  expense  of  more  memory  space,  we  precompute  and  store  the  coefficients 
of  the  vector  field  interpolation  function,  coordinate  transformation  function,  and  the  specialized 
Runge-Kutta  method  during  the  preprocessing  stage. 

As  a  result,  a  cell  record  has  three  coefficient  matrices,  four  node  numbers  and  four  cell 
numbers.  The  four  node  numbers  are  indices  of  nodes  that  comprise  this  tetrahedral  cell.  The 
four  cell  numbers  are  indices  of  cells  that  are  adjacent  to  this  cell.  The  values  stored  in  a  node 
record  include  the  physical  coordinates  of  the  node  as  well  as  the  vector  field  on  the  node.  After 
the  preprocessing  stage,  node  records  become  redundant  and  can  be  deleted  since  the  cell  records 
contain  all  the  information  needed  for  performing  the  particle  tracing. 

Using  our  tracing  method,  each  cell  record  takes 

(3  matrices  +  4  node  indices  +  4  cell  indices ) 

=  (3x(4x3)  +  4  +  4)x4  bytes 
=  176  bytes 

For  a  typical  500,000-cell  data,  about  88  megabytes  are  needed.  If  memory  becomes  a  problem, 
the  matrices  for  the  interpolation  and  the  transformation  functions  can  be  computed  on  the  fly, 
but  the  curl  and  divergence  for  each  cell  then  must  be  stored  at  the  expense  of  much  less  memory 
space,  and  the  memory  requirement  for  each  cell  becomes  96  bytes. 

In  order  to  compare  the  performance  of  the  new  algorithms  with  the  conventional  Runge-Kutta 
methods,  we  have  also  implemented  the  second  and  the  fourth  order  Runge-Kutta  methods. 
Similarly,  to  accelerate  the  tracing  as  much  as  possible,  the  matrices  for  the  interpolation  and 
transformation  functions  are  precomputed  and  stored.  Therefore,  each  cell  record  takes 

(2  matrices  +  4  node  indices  +  4  cell  indices ) 

=  128  bytes 
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However,  the  list  of  node  records  is  needed  during  the  tracing  stage.  On  the  other  hand,  without 
storing  these  two  matrices,  the  memory  requirement  becomes  only  32  bytes  per  cell  record,  and 
24  bytes  per  node  record.  To  cope  with  the  high  memory  requirements  for  visualizing  on  unstruc¬ 
tured  grids,  some  divide-and-conquer  strategies  must  be  taken  to  make  possible  visualization  of 
large  data  sets  such  as  those  with  millions  of  cells. 

5  Test  Results 

To  study  the  performance  of  our  algorithms,  we  compare  experimentally  our  specialized  Runge- 
Kutta  method  (SRK4)  with  both  the  conventional  second  and  fourth-order  Runge-Kutta  methods 
(RK2  and  RK4)  for  integrating  particle  paths.  To  derive  fair  measurements,  as  described  in 
previous  section,  all  the  needed  matrices  are  precomputed  and  stored  for  the  implementation 
of  each  method.  Three  data  sets  are  used  for  our  tests.  The  first  data  set  was  generated 
analytically;  it  contains  68,921  nodes  uniformly  positioned  in  a  cubic  domain,  in  which  there  are 
totally  320,000  tetrahedra.  The  vector  field  on  a  node  is  determined  by  evaluating  three  linear 
functions: 


u-i(x,y,z)  =  — 0.5x  —  6.0y, 

U2(x,y,z)  —  6. Ox  —  0.5  y, 
u3(x,y,z )  =  —  2.0z  +  20.5. 

The  second  data  set  is  the  blunt  fin  data  set  obtained  from  the  National  Aerodynamic  Simula¬ 
tion  Facility  at  the  NASA  Ames  Research  Center.  This  data  set  was  from  a  computational  fluid 
dynamics  simulation  of  air  flow  over  a  fiat  plate  with  a  blunt  fin  rising  from  the  plate  [5].  The 
flow  is  symmetrical  about  a  plane  through  the  center  of  the  fin,  so  only  one  half  of  the  complete 
geometry  is  present.  Note  that  originally  the  computational  grid  was  a  single,  curvilinear,  struc¬ 
tured  block  grid.  We  converted  it  into  an  unstructured  grid  by  splitting  each  hexahedron  into  six 
tetrahedra.  The  resulting  unstructured  grid  contains  224,874  tetrahedral  cells  and  40,960  nodes. 

We  obtained  the  third  data  set  from  the  NASA  Langley  Research  Center.  It  was  from  a 
computational  fluid  dynamics  simulation  of  transonic  flow  about  an  ONERA-M6  wing  with  free- 
stream  Mach- number  0.84  and  3.06  degrees  angle  of  attack  [8].  There  are  287,962  tetrahedral 
cells  and  53,961  nodes  in  this  data  set. 

On  each  data  set,  one  hundred  seed  points  are  randomly  selected.  Then,  streamlines  are 
constructed  by  using  these  seed  points.  The  streamline  constructions  are  stopped  when  either 
the  streamlines  reach  domain  boundaries  or  the  number  of  time  step  exceeds  a  predefined  limit 
(e.g.  1,500). 

Since  the  major  function  evaluations  of  all  the  three  methods  are  of  the  same  kind,  i.e.  matrix- 
vector  multiplication,  we  can  predict  their  performance  by  calculating  the  number  of  function 
evaluations  used  in  these  methods.  For  a  single  step  integration,  only  one  function  evaluation  is 
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Figure  4:  Timing  of  Constructing  Streamlines  on  Data  Set  1 

required  by  using  the  SRK4  method  while  four  function  evaluations  are  needed  if  the  RK4  method 
is  applied  and  two  function  evaluations  are  performed  if  the  RK2  method  is  used.  Theoretically, 
the  SRK4  method  should  be  faster  than  the  RK2  method  by  a  factor  of  2.0,  and  faster  than  the 
RK4  method  by  a  factor  of  4.0. 

The  testing  results  for  the  three  data  sets  are  shown  in  Figure  4,  Figure  5,  and  Figure  6. 
Numbers  are  seconds  and  the  measurements  were  performed  on  a  Sun  SparclO  Model  51  (50MHz). 
Only  the  core  of  the  integration  algorithms  was  measured.  The  test  results  agree  with  our 
analysis;  the  SRK4  method  is  always  the  fastest  method  while  the  RK4  is  the  slowest  one. 

The  average  cost  of  computing  a  single  step  integration  by  using  these  three  methods  are 
listed  in  Table  1.  Note  that  now  the  time  unit  used  is  microsecond.  According  to  the  data 
listed  in  Table  1,  the  speed-up  achieved  by  using  the  SRK4  method  is  slightly  higher  than  2.0 
when  compared  with  the  RK2  method  but  may  be  lower  than  4.0  when  compared  with  the  RK4 
method.  The  lower  speed-up  numbers  and  the  differences  between  different  data  sets  could  be  due 
to  both  the  timing  calculations  and  the  overhead  for  fetching  the  coefficients  of  the  interpolation 
functions,  etc. 

Some  visualization  results  generated  by  using  the  algorithms  described  in  this  paper  are  pre- 
sented  in  Figure  7,  8  and  9.  Figure  7  (a)  shows  streamribbons  visualization  of  the  analytical 
data  set.  From  this  image,  we  can  see  the  streamribbons  spiral  toward  a  critical  point  which  is 
a  saddle  point  in  the  vector  field.  The  streamribbons  are  shaded  using  colors  according  to  the 
velocity  magnitudes.  In  Figure  7  (b),  an  image  of  streamtubes  visualization  is  shown  for  the 
same  data  set  and  using  the  same  initial  seed  points.  This  image  reveals  not  only  rotation  of  the 
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Number  of  Streamlines 


Figure  6:  Timing  of  Constructing  Streamlines  on  Data  Set  3 


exec,  time  (in  fis) 

SRK4 

RK2 

RK4 

Analytical  data 

10.31 

22.20 

32.81 

Blunt  Fin  data 

9.92 

21.33 

34.37 

0NERA-M6  Wing  data 

9.00 

20.90 

32.62 

Table  1:  Execution  Time  of  a  Single  Time  Step 


flow  but  also  expansion  and  contraction  of  the  flow. 

Figure  8  (a)  and  (b)  show  the  streamribbon  and  streamtube  visualization  of  the  blunt  fin  data 
set.  For  both  images,  the  view  is  selected  such  that  the  blunt  fin  is  laid  down  toward  the  viewer 
and  the  plane  surface  becomes  orthogonal  to  the  viewing  direction.  From  these  two  images,  some 
interesting  flow  movements  are  revealed  near  the  leading  edge  of  the  fin  and  the  plane. 

Figure  9  (a)  and  (b)  display  the  streamribbon  and  the  streamtube  visualization  of  the  ONERA- 
M6  wing  data  set.  The  images  show  the  formation  of  a  wing  tip  vortex  caused  by  the  flow 
expanding  around  the  wing  tip  due  to  pressure  differences  between  the  upper  and  lower  surfaces 

of  the  wing. 


The  fourth  order  Runge-Kutta  method  is  the  fundamental  procedure  for  constructing  streamlines. 
A  new  computational  method  has  been  derived  to  speed  up  the  Runge-Kutta  method.  A  closed 
form  formula  is  deduced  to  compute  the  angular  rotation  rate  of  flow  for  making  streamnbbons. 
We  have  also  derived  an  explicit  solution  for  computing  the  radius  of  streamtube  that  is  governed 
by  an  ordinary  differential  equation.  The  performance  of  the  new  methods  were  measured  by 
using  three  different  data  sets  on  a  Sun  SparclO.  The  test  results  match  our  analytic  predictions. 
The  speed-up  currently  achieved  can  be  significant  resulting  in  better  interaction  when  tracing 

a  large  number  of  particles  in  a  large  data  space. 

While  we  have  improved  particle  tracing  calculations,  the  use  of  parallel  processing  can  further 
speed  up  the  tracing  of  a  significantly  large  number  of  particles.  In  addition,  for  data  sets  that  do 
not  fit  into  the  main  memory  of  an  average  workstation,  the  design  of  out-of-core  or  distributed- 
memory  parallel  algorithms  is  needed  just  to  make  visualization  possible.  We  are  currently 
designing  an  out-of-core  particle  tracing  algorithm. 
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